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Abstract 

In this paper we study the Lawrence-Doniach model for layered superconductors, 
for a sample with finite width subjected to a magnetic field parallel to the supercon- 
ducting layers. We provide a rigorous analysis of the energy minimizers in the limit 
as the coupling between adjacent superconducting layers tends to zero. We identify 
a unique global minimizer of the Gibbs free energy in this regime ("vortex planes"), 
and reveal a sequence of first-order phase transitions by which Josephson vortices 
are nucleated via the boundary. The small coupling limit is studied via degenerate 
perturbation theory based on a Lyapunov-Schmidt decomposition which reduces the 
Lawrence-Doniach system to a finite-dimensional variational problem. Finally, a lower 
bound on the radius of validity of the perturbation expansion (in terms of various 
parameters appearing in the model) is obtained. 
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1 Introduction 



In 1971 Lawrence and Doniach [[LaDo 71| introduced a Ginzburg-Landau type model for 
superconducting materials with a planar layered structure. In this model, the superconduc- 
tor occupies an array of parallel sheets with insulating material acting as a buffer between 
the sheets. While this model was originally proposed to study layered structures artificially 
produced by successively deposing thin planar sheets of superconducting metal with organic 
separators, it has received renewed attention due to the discovery of high temperature su- 
perconductors. Indeed, nearly all of these high-T c materials are crystals with a distinctly 
layered structure, consisting of copper oxide superconducting planes stacked with intervening 
insulating (or weakly superconducting) planes. 

In this paper we will consider the case of a layered superconductor in a uniform magnetic 
field imposed parallel to the superconducting planes. We assume that there are a finite 
number of superconducting sheets, each parallel to the xy-p\a.ne, with uniform separation p. 
We assume that the external magnetic field is applied along the y-direction, H = Hy. We 
will take the planes to be of infinite extent in the y-direction, and assume the local magnetic 
field will be everywhere independent of y and point in the y-direction, 

K(x,y,z) = h(x,z)y. 

The vector potential A may then be chosen to lie in the xz-plane, 



A(x, y, z) = A x (x, z)x + A z (x, z) z, 



h = curl A 



' dA T dA, 



dz 



dx 



We assume that the sample has fixed width 2L in the x-direction, and hence the supercon- 
ducting sheets are described by the stack of parallel planar strips 

S^v : —L < x < L, — oo < y < oo, z = z n := np, n = 0, 1, . . . , N. 



Since each sheet XV is superconducting it carries a (complex- valued) order parameter 
i/j n {%), n — 0,...,N. We choose units in such a way that = 1 represents a purely 
superconducting state. The Lawrence-Doniach model is then formulated in terms of the 
following Gibbs free energy functional: 

\ 2 



g r (ifj nj A) 




— -iA x \ ip n 
. dx / 



tp n - ip n -.i exp ( i I A z (x, s) ds 

%n— 1 i 



dx 



dx 



+ 



L fN P (QA dA 2 



dz dx 



H dzdx 
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where r is the interlay er coupling parameter (or Josephson coupling parameter.) We have 
chosen units such that the in-plane penetration depth X ab — 1, k — \ a b/£,ab is the Ginzburg- 
Landau parameter, (A a &, £ a {, are the in-plane penetration depth and correlation length, re- 
spectively,) and the magnetic fields are measured in units of H c / k, where H c is the thermo- 
dynamic critical field. (See [ [T 96|| .) 



The coupling between the superconducting planes given by the second sum in Q r simulates 
the Josephson effect, by which superconducting electrons travel from one superconducting 
sheet to another by quantum mechanical tunnelling. We will see this explicitly in the Euler- 
Lagrange equations, where the currents in the gaps between planes will be determined by 
the sine of the gauge-invariant phase difference. The interlayer coupling parameter r gives 
the strength of the Josephson effect. In our units, 

2 



XjK 2 p 2 



where A j is the (non-dimensional) Josephson penetration depth. In an anisotropic Ginzburg- 
Landau model, Aj gives the effective mass ratio which determines the degree of anisotropy. 
For highly anisotropic superconductors Aj is very large, and hence when up ~ 1 we may 
treat small parameter. 

Due to the layered structure one expects these materials to be highly anisotropic. A first 
attempt to model layered superconductors is by an anisotropic Ginzburg-Landau model, 
which treats the sample as a three-dimensional solid with anisotropic material parameters. 
For certain materials and temperatures close to the critical temperature T c this approxima- 
tion seems valid, but for the most anisotropic superconductors the anisotropic Ginzburg- 
Landau model does not give a good qualitative or quantitative description of experimen- 
tal observations. For example, when the sample is subjected to a sufficiently strong mag- 
netic field oriented parallel to the superconducting planes Kes, Aarts, Vinokur, and van der 



Beek ||KAVB 90| observe a transition betwen "three-dimensional" behavior (governed by the 
anisotropic Ginzburg-Landau model) and "two-dimensional" behavior at a critical tempera- 
ture Tco below Tc- In the two-dimensional regime the superconducting planes decouple and 
the applied magnetic field penetrates completely between the planes, virtually unscreened 
by the superconductor. Despite the penetration of the field, superconductivity within the 
planes is not destroyed even in very strong applied fields. This "magnetically transparent" 
state is inconsistent with the anisotropic Ginzburg-Landau model, where magnetic fields of 
moderate intensity are largely expelled from the bulk except for an array of isolated vortices 
(the "Abrikosov lattice".) In addition the Ginzburg-Landau model predicts the breakdown 
of superconductivity when the applied field penetrates the material completely, which occurs 
when the field exceeds a critical value H c2 . 

We note that Chapman, Du, & Gunzburger ||CDG 95f have proven that solutions of 
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the Lawrence-Doniach model converge to solutions of the anisotropic Ginzburg-Landau 
model (and in particular the convergence of energy minimizers) under the limit p —>■ with 
k, Xj fixed. This limit does not correspond to our "two-dimensional" regime, since it would 
send r — > oo, corresponding to a strong coupling between adjacent superconducting layers. 
Indeed, we observe that the non-dimensional separation distance p of our model is related 
to the (dimensionally dependent) physical separation p via p = p/X a b- Since \ a b(T) — > oo 
as T — > T c , the limit T — > T c effectively corresponds to p — > (with n, Xj,p fixed) in our 
units, and therefore the Chapman, Du, & Gunzburger limit can be interpreted as letting 
T — ► T c . This is consistent with the observed "dimensional crossover" to the anisotropic 
Ginzburg-Landau regime for temperatures near T c . 

In this paper we will study the minimizers (and low-energy solutions) of the Lawrence- 
Doniach system for r near zero and thereby analyse the structure of the resulting "transparent 
state." The crucial observation is that when r = the planes decouple, and the energy may 
be minimized explicitly by solving simple first order equations. Even after gauge symmetries 
have been removed the r = problem exhibits an additional symmetry, corresponding to an 
iV-dimensional torus action (where N + 1 denotes the number of superconducting planes,) 
and thus the minimization problem at r = degenerates on a finite dimensional manifold 
in function space. We can think of the r = problem in analogy with the self-dual point of 
the Ginzburg-Landau model (k = l/y/2), where minimizers satisfy a first-order Bogomolnyi 
system (in addition to the usual second order Ginzburg-Landau equations,) and the same 
minimum energy is obtained by any configuration of vortices. 

When r ^ this symmetry is broken and a Lyapunov-Schmidt decomposition reduces 
the problem of finding solutions with r ~ to a finite dimensional variational problem on 
the degenerate manifold. The minimum value of energy is 0(r), and we indeed recover the 
"transparent state" observed in experiments. The local magnetic field h(x,z) = H + 0(r) 
inside the sample, and superconductivity is hardly affected in each plane, \if) n {%) I = l — 0(r). 
In particular the order parameters are never zero: "vortices" correspond to local maxima 
of the local magnetic field, and lie between the layers. In the physics literature these are 
referred to as Josephson vortices, as opposed to the Abrikosov vortices typically observed in 
the Ginzburg-Landau theory. 

The finite dimensional reduced problem may be solved explicitly to determine the exact 
geometry of the Josephson vortex lattice. For a finite sample in x, z the minima of energy 
form "vortex plane" configurations, in which the local magnetic field is uniform in z. The 
vortices are not separated, but line up vertically at the local maxima of h(x, z) = h(x). (See 
figure 1.) As the external field H is increased vortices are nucleated at the edges, by a first- 
order phase transition. Formal asymptotic expansions for these solutions were calculated by 
Theorodakis ||Th 90|| , and Kuplevakhsky [|K 99|| claimed that they were the only solutions of 
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the Lawrence-Doniach system. On the contrary, we find that there are exactly 2 N solutions 
with energy 0(r): two represent vortex planes (one stable and the other unstable), and the 
others (unstable) lattices of various geometries. (See Theorem |4.1| .) 

Again, we note the distinction with the Ginzburg-Landau model: the geometry of the 
Abrikosov lattice was determined by numerical comparison of a finite number of possible 
lattice geometries. For the Lawrence-Doniach model in the small coupling limit we are 
able to identify the absolute minimizer (and all low-energy solutions of the Euler-Lagrange 
equations) explicitly and rigorously. This is a direct benefit of the discrete nature of the 
model. 

The basic idea that an infinite dimensional variational problem is actually governed by 
a finite dimensional one in some parameter limit is a common one in analysis. Indeed, it 
appears in such diverse contexts as the location of Ginzburg-Landau vortices as k — > oo 
(Bethuel, Brezis, & Helein [ [BBH 94f| , Bethuel & Riviere ||BR 95|| ,) spike-layer solutions (for 
example, Li & Nirenberg |[LN 98|| , Gui ||Gu 96|| , or Wei ||W 98|| ), and blow-up for critical expo- 
nent problems (for example in Bahri, Li, & Rey | [BLR 95| , or Rey |[Re 91|| .) These examples 
are of singular perturbation problems, though. The transparent state arises as a degenerate 
regular perturbation of the r = problem, and hence it is more closely related to the work 
of Ambrosetti, Coti-Zelati, & Ekeland ||AmCzE 87|1 and Ambrosetti & Badiale [[AmBa 98 



on 



homoclinic solutions of Hamiltonian systems and the Poincare-Melnikov functions. 

Analytically our results are unambiguous: for any choice of the other parameters (L, k, 
N, H, p) we can choose r sufficiently small so that the vortex planes configuration minimizes 
the free energy. However, in a real superconductor r is not infinitessimally small, which raises 
the question of how the interval of validity of the r-expansion is affected by the values of 
the other parameters in the problem. For example, in the order r term in the expansion 
of the solution (see ([34]), for example) we observe " secular" terms, that is factors which 
become large without bound as the length of the interval increases, and in general remark 
that the coefficients increase with L, k and decrease with H. We address this question in 
section 5, where we produce a lower bound for the interval of validity as a function of N, 
L, k,, and H. We discover that this interval is independent of N and indeed increases with 
decreasing L, k and increasing H. For the high-T c superconductors, k is large and typical 
macroscopic sample widths L are very large compared with the in-plane penetration depth 
X a b. This indicates that our analysis may be more applicable to experiments with highly 
anisotropic organic or synthetic multi-layer superconductors, where the material parameters 
are significantly different from the high-T c crystals. 

In a subsequent paper ||AlBeBr "OOf we address the question of minimizers of the Lawrence- 
Doniach energy in very large samples, by considering periodic solutions in an infinitely wide 
sample. In fact, in the periodic case we find a different geometry for the energy minimizing 
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configuration! The periodic solution with the least energy is a lattice with period two in n, 
forming a diamond pattern of Josephson vortices, proposed by Bulaevskii Sz Clem [[BuCm 91|j . 
In the periodic case, the role of L is now taken by the period, which for minimizers will 
decrease as the applied field H increases. Therefore the small r expansion will have a 
large range of validity in sufficiently large fields H, and the result should better describe 
experiments with high-T c superconductors. 

Acknowledgement. The authors thank A. Bahri for his helpful suggestions for determining 
the region of validity of the perturbation methods. SA and LB wish to thank R. V. Kohn 
and P. A. Deift for inviting them to spend a semester at the Courant Institute in Fall 1999, 
and for the interest and encouragement they have shown throughout the years. 

1.1 Variational setting 

We begin with the following basic energy estimate, which legitimizes the simplification 
|^ n (a;)| ~ 1 for small r. Indeed, in the physics literature the approximation \tp n \ ~ 1 is 
widely assumed to hold: see for example Bulaevskii ||Bu 73|| , Clem & Coffey [CmCo 90|| , or 



Bulaevskii & Clem PuCm 9"1|. In particular, we see immediately that in the regime r ~ 



there will be no "Abrikosov vortices" in the sense that the order parameter is never zero. 
Proposition 1.1 Let D := [—L,L] x [0, Np\. For any r > we have 

mf{^-g r (ip n ,A): ip n eH\[-L,L]), A e H\D; R 2 )} < 2Np(L + -\-)r. 
H£ pn 

Moreover, there exist constants ro > and C = C(N, L, k, H,p) > such that for all r G 
[0, ro) the minimum is attained by (ip n ,A) with \ijj n (x)\ > 1 — CV 1//2 > for all x E [—L,L]. 

Proof: We choose a test configuration, ip n (x) = exp{inpHx}, A = (Hz,0) to estimate the 
free energy, 



4?r N r L / l \ 

—QM W A)=rpJ2 (1 - cos{pHx)) dx < 2N I L H — — pr. 



The fact that the minimum is attained can be easily proven once an appropriate choice of 
gauge has been made: see ||CDC 95f for details. 



It remains to show that the order parameters lie near the unit circle. We recall the 



"Diamagnetic inequality" (see for example p. 174 of | LL 97|1 ), which states 



| (x) < |(V — iA)f\(x) (almost everywhere) 
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for every / G L 2 oc with (V — iA)f G Lf oc . Using the elementary inequality (1 — \ip n \) 2 < 
(1 — |^ n | 2 ) 2 ) and the energy bound we obtain 



N N L ( 1 

Elli-Wll^i < E/ (i-W 2 ) 2 + ^ 

n=0 n=0 J ~ L \ K 



dx 



N r L ( Id \ 

< E / l(i-M 2 ) 2 + - 2 \(^-iAM 2 ) dx 



n=0 



4:71 I 1 

< jpg r {i/> n ,A)<2N\L + — )i>r. 



k 2 dx 
Hp, 



Here we choose an equivalent norm for H l ([—L, L\), 



By the Sobolev embedding we have for each n = 0,...,N that |^ n (x)| > 1 — CV 1 / 2 with 
constant C depending on N, L,p, k, H . 



The above proposition suggests the use of polar coordinates for ip n in order to deal more 
directly with the phase of the order parameter, which plays the essential role in Josephson 
coupling. We define f n ,4> n via i/) n ( x ) — fn( x ) exp(i0 n (x)), and note that <p n is well-defined 
only up to an additive integer multiple of 2n. We then define our new free energy in terms 
of the variables (f n ,(p n ,A) to coincide with jjkG r , that is: 



n r (f n , (f) n ,A) = P E / -Sfn - l) 2 + -Unf + ~M - Zn)) 2 f 2 
n=0 J ~ L VA K K 

r N f L 

+ nP E / (fn + fn-1 ~ Vnfn-1 COs($ n , n _i)J dx 



dx 



+ 



1 ,L f Np ( QA dA 7 



K J-L JO 



dz dx 



— H ) dz dx, 



where 



$ n ,n-l(x) ■= 0n - 4>n-l ~ / A z (x, z) dz, 



Zn-l 



is the gauge-invariant phase difference, z n = np and r 



X 2 jK 2 p 2 



. When f n — \ip n \ is bounded 



away from zero the condition ip n G H 1 ([—L, L]) is equivalent to both /„, <p n G H 1 ([—L, L]). 
We first define a base space, in which Q r will be a smooth functional: 



£ :-- 



{f n) <f> n , A): f n G H\[-L, L]), <t> n G H\[-L, L}), n = 0, . . . , N 
A=(A X ,A Z )EH\D,R 2 ) 
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We remark that we should really work in the convex subset of 6 with f n > for all n = 
0, . . . , N, but Lemma |1.1| already guarantees that the solutions we will find will have f n ~ 1. 
Furthermore, Proposition 4.5 of [|CDG 95|| asserts that f n (x) = \ip n (x)\ < 1 for any solution 
of the Lawrence-Doniach equations. 

We note that Q r is a smooth (C°°) functional on S, and that variation with respect to 
each of its arguments gives the Lawrence-Doniach system. First, we denote by 



Vn 



~ Ar (^5 z n)), 



X 



the supercurrent velocity. Then, variation of Q r with respect to /„ (for each n = 0, . . . , N) 
yields 



f I (fn- 

COS $n,n— 1 fn+1 COS $n+l,n 

2f n ), n^0,N; 



K 



| (fl cos $i,o - fo) , 
t | (fN-1 COS $AT iA r_i - /at) 



with boundary condition f^(±L) = 0, n = 0, . . . , N. 

Variation with respect to A x produces the following equation in weak form: 



n = 0; 
n = N, 



(2) 



dh 

dz 



(x, z) = for all (x, z) G (— L, L) x (z n -iz n ), n = 0, . . . , N, 



h(x, z n +) - h(x, z n -) = -pfl{x){(f)' n - A x {x, Zn)), n = 1, . . . , iV - 1, 
h(x, z N -) = pf N (x){(j)' N - A x (x, z N ), h(x, 0+) = -p/ o 2 (x)(0' o - A x (x, 0). 



(3) 
(4) 



In other words, h is independent of z away from the superconducting planes, and supercur- 
rents in the SC planes create jump discontinuities. The dependence of h(x, z) on x in each 
gap is determined by (compactly supported) variations of A z : 



= -^-fn(x)f n -i(x) sin (a;), 



if — L < x < L and 

z n -i < z < z n ,n = 1, 



(5) 



with boundary condition h(±L, z) = H . It is therefore natural to define 
h(x, z) = h^ n \x), when — L < x < L, z n -\ < z < z n , n = 1 



,N; 



,...,N, 



(6) 



with determined by the ordinary differential equations (|5|) together with the boundary 
condition h {n \±L) = H. 
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(7) 



Variation with respect to <p n produces the following current-conservation laws: 

' f L/n/n-i sin $ n ,„_i - /„ + i/„sin$„ + i >n ], n = 1, . . . , N - 1; 

-|/i/osin$i j0 , n = 0; 

k |/jv/jv-i sin $ N>N -i, n = N, 

with boundary condition 

/n 2 (^)(0n(^)-^(^^))=O, x = ±L, n = 0,...,iV, 

which expresses the physical fact that current should not flow past the edge of the material. 
Note that (|?p does not yield any new information, since it can be obtained by differentiat- 
ing the jump condition ([|) and substituting from (||). This is not surprising, since gauge 
invariance implies a nontrivial relationship between the <p n and A. Indeed, denoting the 
supercurrents in the planes by 

■= V n fl (8) 

and the Josephson current in between the n — 1 and the nth planes by 

£ n) := r -^ 2 vfnfn-x sin $ n , n _i, (9) 

equation (0) is a semi-discrete version of the classical continuity equation divj = 0, and 
gives the conservation law corresponding to the U(l) gauge invariance in accordance with 
Noether's Theorem. 

We are most interested in the gauge invariant "observable" quantities, which enter di- 
rectly into the free energy: the density of superconducting eletrons f n , the supercurrent 
velocity V n , the gauge-invariant phase difference $ n>n _i, and the local magnetic field h(x, z). 
A very useful formula for <& n ^ n _i(x) can be obtained by applying Stokes' Theorem in the 
rectangle R = (0, x) x (z n -i,z n ): 

Vif^AK-y^ + pf^^^ + VifO) n = l,...,N. (10) 
Jo Jo 

From these equations we easily verify the smoothness of observables associated to weak 
solutions of the Lawrence-Doniach system: 

Proposition 1.2 Suppose (f n ,<p n ,A) G £ are critical points ofQ r . Then f n ,Vn,^n,n-i £ 
C°°([-L, L\), and h G C°°([—L, L\ x (z n . l} z n )), n = l,...,N. 

Of course, the regularity of the non-gauge-invariant quantities <p n and A depends on the 
choice of gauge. 
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Remark 1.3 Some authors have (correctly) pointed out that it is not physically correct 
to impose the external field via a Dirichlet condition h = H on the boundary dD of the 
sample. A more appropriate model for the effect of an external field is obtained by placing 
the superconductor D in a larger region D D D, (with D = R 2 possibly) and including the 
field energy in D \ D in the calculations of the free energy, 

Qr(fn, 0n, A) : = fi P (/ n , <p n , A) + ff (curl A — Hf dxdz. 

JJd\d 

For example, our sample D may be lying at the center of a long cylindrical solenoid of large 
radius (whose interior cross-section is a large disk D.) All of the preceding analysis can then 
be carried through for Q r by choosing an appropriate space (see Rubinstein & Schatzman 
RuS 98|| for a discussion of the correct setting when D = R 2 ,) but it is easy to verify that 



the Euler-Lagrange equations yield h(x, z) = H for (x, z) G D \ D. This is due to the 
two-dimensional ansatz: V x (g(x, z)y) = in a domain implies that g is constant there. 
Consequently we may use the simpler form of the energy Q r with no loss of generality or of 
physical relevance. 

There is a large degree of degeneracy of Q r in £ due to the gauge invariance: if x G 
# 2 ([-L,L] x [0,Np])), and 

fn = fn, (f) n (x) = (p n (x) - x(x,Z n ), A = A - V%, 

then Q r (f n , <fi n , A) = Q r (f n ,(f) n ,A). As usual, we eliminate this troublesome degeneracy by 
fixing a gauge. The most convenient choice is the Coulomb gauge, which allows us to control 
the H l norm of the vector potential A by its curl. We define a subspace of £ to incorporate 
this choice of gauge, 

£ g := {(/„, 4> n , A) G £ : / 4>o(x) dx = 0, div A = in D, and A ■ n = on dD.} 

J — L 

This choice is made with no loss of generality: 

Lemma 1.4 For every (f n ,<f> n ,A) G £, there exists x G H 2 (D) so that (f n ,(p n — x(', z n),A — 

V X ) e £ g . 
We also have: 

Lemma 1.5 For every (f n ,(j) n ,A) G £ g 

1(D) ^ M)||^LU1^L|| L 2(.D), 



H HD) < C ||curlA||| 2(jD) , (11) 



where 

4 L 2 N 2 p 2 



Co 



1 + 



7T 2 N 2 p 2 + AL 2 



'12) 
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Proof of Lemma \1.J\ : Given (/„, n , A) G S we may solve the linear Neumann problem, 

= div A for x G D, 
A ■ n for x G dD. 

By the divergence theorem there exists a unique (up to constants) weak solution, x £ H X (D). 
Since A G i/ 1 we have x £ Hf oc (D) by standard regularity theory. Since the domain 
is polygonal we must be more careful to determine the regularity at the corners, but by 
Theorem 1.5.2.4 of Grisvard ||Gr 85|] we may find an H 2 (D) function whose normal derivative 



coincides with A ■ n on dD. Then Theorem 4.3.1.4 of [|Gr 85|| provides global regularity, 
X G H 2 (D). By subtracting a constant from x such that 4>o(x) — x(%, 0) has average zero we 
achieve the desired gauge change. 



Proof of Lemma \1.3j : Assume (f n ,4> n ,A) G £ g . We now solve the Dirichlet problem, 

J Ar] = curlyl for x G D, 

1 rj = for x G dD. [ ' 



By Theorem 4.3.1.4 of ||Gr 85|| the unique solution rj G H 2 (D), and A = (d z rj, —d x rj) G 
H\D;R 2 ) with div (A— A) = 0, curl (A — A) = 0, and (A- A) -n = on dD. Hence A — A, 
and the explicit solution of the Dirichlet problem in a rectangle (via Fourier analysis) provides 
the constant Co > such that: 



m(D) < \\v\\h*(d) < c o\\^v\\h(D) = C ||curll|| 





2 Minimization at r = 0. 

When r = 0, the superconducting planes decouple, and we may solve the minimization 
problem explicitly. The solution of the problem r ~ will require some detailed second- 
order information on the minimizers at r = 0, so we first establish a functional analytic 
setting for the equations. We exploit the Hilbert manifold structure of £ g , and regard the 
first variation of energy as elements of the tangent Hilbert space 

E = TS g = |(u n ,u n ,a): u n , v n G L]), J ^ v dx = 0, 

aeH 1 ^), diva = in D, a ■ n = on dD.} , 



S. Alama, A.J. Berlinsky, & L. Bronsard 



12 



and the second variation as a self-adjoint operator on E. We first introduce an equivalent 
inner product: for (u n , v n , a), (U n , V n , A) G E, let 



),(U n ,V n ,A) 



J ^2u n U n -\ — -[u' n U' n + v' n V^] + fnKij cfoH — 2 JJ cur ^ ' curlAc/rcck. 

Then we define the gradient of fl r (f n ,(f) n , A) via 

Vft r (/ n ,0 n ,A),( )) = DQ r (f n , 4> n , A)[u 



Proposition 2.1 VFrf/i r = 0, 

(a) inf {Q (f n , 4> n , A) : (/„,</>„, A) G £ 9 } = 0. T/ie minimum value is attained, and the set 
of all minimizers coincides with the hyperplane S, defined by the set of (f n , <p n , A) G S g 
such that: f n = 1; 

rx 

<p n (x) = a n + AJs,z n )ds, a = 0, a n G R, n = 1, . . . , N; (14) 
J o 

A = (r^, — 773;) , where Ar] = H in D, and T]\go = 0. (15) 

In particular, the gauge invariant phase difference is given by $° n _ 1 = 5 n + Hpx, where 
5 n := a n - a n -i, n = l,...,N. 

(b) For any element s = (f®,^, A°) G S, the linearized operator D 2 Qq(s) : E — > E 
defines a Fredholm operator with index zero. Moreover, 

T S S = ker D 2 tt (s) ~ R N . (16) 

Note that T := T S S is independent of s G 5. Note that 5 may be parametrized by either 
(ai, . . . , ttTv) or (^i, • • • , <V) ^ R-^ — We abuse notation and write s = s(cti, . . . , a^v) or 
s = s(6 n , . . . , 5jv) to emphasize the dependence. 

Proof: Being a sum of non- negative terms we clearly have inf,f 9 Q (f n , <p n , A) > 0. The 
infimum of zero will be attained at (/°, 0°, A ) if and only if they solve the following first- 
order equations in our space S g : 

f n {x) = l, <P' n (x)-A x (x,z n )=0, cm\A = H. (17) 

Note that by Lemma |1.4j the last equation is uniquely solved for A G H 1 (D), with solution 



as in (|i~5P . By the trace theorem A x (x,z n ) G H l l 2 {\— L, L]), and therefore, <p n is uniquely 
determined by integration, except for the (N + 1) constants of integration ot n . Note that 
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r) is even in x, and hence Jq A®(t, z n ) dt is an odd function. In particular, a n gives the 
average value of </>° and thus a = is fixed by the definition of the space £ g , leaving N free 
parameters («i, . . . , «tv) to parametrize the solution set. The explicit form for $° n _ l is then 
obtained from fllPl) . This completes the proof of (a). 

Writing the linearized operator around a solution 0°, A°) as a quadratic form, 



N r L r 1 1 1 

D 2 n Ul<PlA»)[u n ,v ni d] = pY / J_ L [2ul + - 2 [u' n } 2 + - 2 K-a x (x,z n )} 2 j dx 



+ ~ r // I cur l ^ I 2 dx dz 



K 2 JJD 



{(u n ,v n ,a), (u n ,v n ,a)) 



N f L ( 2 ] 
PJ2 / ^(^(a;,^)!^) (18) 



where (m„, v n , a) G 

First note that by the second identity in flTB|), D 2 Q (f®, 0°, A ) is of the form identity plus 
compact (since the trace embedding of a — >• o(-, z n ) is compact from H 1 ^) to L 2 ([— L, L]).) 
Next assume that (u n , v n , a) G ker D 2 Q (f®, 0°, A ). Clearly u n = for all n = 0, . . . , AT. By 
Lemma |1.5| we must have a = 0, and hence t> n are constant. Since the definition of the space 
forces t>o(x) with mean zero, we are left with iV free parameters, and 

kerD 2 fi o (/°,0°,^) = { a) G E : (u n ) n = 0, a = 0, 

(u , ■ ■ ■ , v N ) = (0, Ci, . . . , cat), ci, . . . , G R} 

= T S S. 





3 Degenerate perturbation theory. 

We now perturb away from the degenerate minima of Qq, using a variational Lyapunov- 
Schmidt procedure. This method has been used by Ambrosetti, Coti-Zelati, & Ekeland 
AmCzE 87| , Abrosetti & Badiale |AmBa 98| , Li & Nirenberg | |LJN 98| | (and many others) 



in a variety of situations involving heteroclinic solutions of Hamiltonian systems and in the 
semiclassical limit of the nonlinear Schrodinger equation. 

We now proceed according to the usual Lyapunov-Schmidt reduction. Since S is a 
hyperplane, T = T S S is independent of s G S. Let W = T- 1 , so any (f n , <p n , A) G S g admits 
the unique decomposition (f n , <fi n , A) = s + w with s G S, w G W, and any U := (u n , v n , a) G 
E decomposes uniquely as U = t + w with t G T, w G W. We denote the orthogonal 
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projection maps P : E — ► T, P -1 : P — > W so that P/7 = t, P L U = w whenever U = t + w. 
As in the previous section we interpret the first variation Vfl r (/„, <f> n , A) as an element of E 
itself, and project the equation Vfl r (f n , (fi n , A) = into the two linear subspaces T and W, 

Pi(r,s,u>) := P [VQ r (s + w)} = 0; (19) 
P 2 (r,s,w) := P x [VO-(s + ^)] = 0. (20) 

By Proposition [TT] (b), the second equation can be solved uniquely for w = w(r,s) in a 
neighborhood of S for r small, using the Implicit Function Theorem. Because our functional 
Q r is smooth we can expand w(r, s) in powers of r. Note that f2 r (s(«i, . . . , ocn) + w) is 
periodic in (ax, . . . , a at) so that we may think of 5 as a (compact) iV-torus. Therefore the 
expansion will be uniform in s. We obtain the following variant of Lemma 2 of ||AmBa 9"8|| : 

Lemma 3.1 There exist constants r > and 5 > 0, depending on N,L,k, and H, and a 
smooth function 

w = w(r, s) : (— r , r ) x S — > C E 

such that: 

(i) There exists smooth functions W\ 1 W2 such that 

w(r, s) = rwi(s) + r 2 w 2 (r, s) 
for all |r| < r and for all s G S; 

(ii) P ± [Vn r (s + w(r,s))] = 0. 

(m) Conversely, if P ± [VQ r (s + w)] = for some r G (— r ,r ) and it) G W rai/i < 5, 

then w = w(r, s). 

(iv) For any choice of Lq, Kq, Hq > the constant ro may be chosen uniformly for all N > 1, 
1 < L < L , 1 < k, < k,q, and H > H . 

Parts (i)-(iii) follow easily from the Implicit Function Theorem. The dependences on the 
various parameters is more delicate: we provide the full proof in section 2.5. 

We define 

S r := {s + w(r, s) : s G S} C £ g . 

S r is a smooth manifold smoothly diffeomorphic to the hyperplane S. The important role 
played by S r is that it is a natural constraint for Q r (see Lemma 4 of ||AmBa 9"8fl ), and hence 
the equation fll9|) may be solved variationally: 
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Lemma 3.2 

(a) If (f n , (f) n , A) e S r satisfies D(tt rlSr )(f n , (f> n , A) = 0, then Vfi r (/„, (f> n , A) = in E. 

(b) There exists e = e (N, L, k, H) > such that if (f n , <fi n , A) G £ g is a critical point of Q r 
with 

tt r (fn, <frn, A) < € , (21) 

then (f n ,(f) n ,A) e S r . 

(c) There exists fo = ro(N, L, K , H) such that for all < r < fo, inf,f 9 Q r = infs r ^r- 

Proof: The assertion (a) is exactly Lemma 4 of ||AmBa 9"5| . To prove (b), assume that 
{fni'PniA) e £ 9 is a critical point of f2 r with energy bounded by (pll) . As remarked in the 
beginning of the section, (/ n , n , A) decomposes uniquely as (f n , (fi n , A) = s + w, s = (/° = 
1, 0°, A ) 6 5, w 6 VF. Now we use the energy bound to estimate the distance between 



(f n , (p n , A) and s: by Lemma |L5 



1^4 — AoIIh 1 < C 1| curl A — curly4 || L 2 



n=0 - 



■(/;) 2 + (i-/n) 2 

•£ 1 



C 1| curl A — H\\ L 2 < C e K 2 ; 
N t-L r 1 



if'n) 



f 

J n 



dx < eo; 



dx < p E / 

n=0 ^~ 

" f 1 * 1 , N f L 1 

PZ / — = PH — (4>' n ~ A x0 (x,Z n )) 2 dx 

N fL / 1 

+— (A^x,^) - A x0 (a;,2; n )) 2 ^ dx 



< 2e Q + C 1 \\A-A \\ m <C 2 e , 



where Ci comes from applying the Trace theorem (see Lemma p.2| in section 5) in the last 
line. Finally, since [(f n , (fin, A) — s] 6 W, each (fi n — (fi n has mean value zero. Therefore, 
the ff 1 ([— L, L])-norm of the difference is controlled by the difference of the derivatives, as 
estimated above, and we may choose eo small enough such that 



dist ((/„, (fi n , A), (/no, (fino, v4 )) < C6 < 5, 



(22) 



where 5 = 5(N, L, k, H) is given by Lemma |3.1] . It then follows by assertion (iii) of Lemma |3]T 



that if (f n , (fin, A) is a critical point of the Lawrence-Doniach system with the given energy 
bound it must lie on S r . This completes the proof of (b). 

To prove (c) we note that Proposition |1.1| implies that inf^ g Q r < 2Np(L + jj-)r, and 
hence we can choose f < r such that 2Np(L + ^)r < e . 
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Remark 3.3 We note that we cannot make the same statement about the N, L, k, H depen- 



dences of f in (c) of Proposition |3.2| as we make for r in Lemma |3.1| . The uniform bounds 
on r are possible because of local estimates on the solution set of fl20|) and a continuity 
argument from r = 0. (See section 5.) We have no such control on the distance 5 from 
S to solutions which do not lie on on the manifold S r . This would entail uniform (in the 
parameters) global (i.e., non-perturbative) information about the solution set, which energy 
bounds do not provide. This leaves open the possibility of an interval r G (ro,ro) for which 
the solutions on S r continue to exist and are represented by a perturbation expansion in r, 
but the absolute minimizer might not be an element of this family. 
<> 

In conclusion, we have achieved a complete finite-dimensional reduction of our problem, 
for small r. That is to say, when < r < f all low energy solutions of the Lawrence-Doniach 
system can be found on the N- dimensional surface S r . Moreover, an explicit form for these 
solutions may be determined by a simple procedure of Taylor expansion of the equations and 
energy in powers of r, as is legitimized by Lemma |3.1| . 



4 Vortex planes. 

We now apply the theory of the previous section to determine the minimizer (and other 
stationary states) of the Lawrence-Doniach energy, for r << 1. 
We summarize our results in the following: 

Theorem 4.1 Assume that sin(HpL) ^ 0. There exists T\ = ri(N,L,K,H) such that for 
every r with < r < r±, the global minimum of free energy is attained by the vortex plane 
solutions, given asymptotically by ft23j)-ffij\) below. Moreover, when < r < r 1 Q r admits 
exactly 2 N physically distinct critical points with energy bounded as in j\2l{). 

4.1 Minimizing Q r \s r 

By Proposition ^| we seek critical points of the finite dimensional variational problem Q r \s r - 
Using Lemma [Tl] (i) we expand a point (f n , n , A) e S r as 

(/„, n , A) = s + w(r, s) = (/°, 0°, A ) + r(u n , v n , a) + 0(r 2 ), 

where (/°,0°, A ) G S solve the first-order system (p]), and the error term is uniform over 
s G S. We observe that fl r has the form, 

tt r {f n ,<f>n,A) = n (f n ,(f) n ,A)+rT(f n ,(j) n ,A), with 
r(/ n ,0 n ,i) = (/„ 2 + /„ 2 -i-2/ n /„- 1 cos($ n , ri _ 1 )) dx 

1 71=1 J ~ L 
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Since Q r (s + w(r, s)) is a smooth function of s = s(5i, . . . , Sn) and r, it admits a Taylor 
expansion of the form, 



d 

Q r (s + w(r,s)) = Q (s) + r —QJs + w(r, s)) 

or 



+ 0(r< 

r=0 



r r(s) + vn (s) 



dw , 



+ 0(r 2 ) = rr(s) + 0(r 2 



rp^2 [I — cos(S n + Hpx)]dx + 0(^ 

n=l 



= 2 P (iVL - Sin( ^ L) f] cos r + 0(r 2 ), (23) 

with remainder term uniform for s G S. 
Define G:RxR w ^Rby: 

G(r, 5i, • • • , S N ) := tt r (s + w(r, s))/r, s = s(5 u 5 N ). 

Then G is a smooth function which is periodic in each coordinate S n , and 

G(0, <fo) = 2iV>L - 2sin( ^ L) £ cos 5 n . 

^ n=l 

Note that when 

sm(HpL) 

G(0, s) possesses exactly 2 N critical points for (Si, . . . , S n ) G JC :— H N / (2nZ) N corresponding 
to any permutation of 

S n G {0, 7r} mod 27T. 

It is easy to see that each is a non- degenerate critical point of G(0,s). By the Implicit 
Function Theorem, there exists f i > such that for all r G (0, fi), and for any critical 
point (<5*, . . . , 5^) of G(0, Si, ... , Sn) there exists a unique critical point (<5i(r), . . . , Sn(t)) of 
G(r, <5i, . . . , Sn), with 

(5 1 (r),...,5 N (r)) = (S* 1 ,...,S N ) + 0(r). 

Since G is periodic in each S n we may also conclude (via a compactness argument) that these 
are the only critical points of G(r, Si, ... , Sn) for small r. Since f2 r |s r = rG(r, Si, ... , Sn), by 
Lemma [3.2| , therefore, for all r with < r < min{fQ,fi} := r\ and whenever sin(HpL) ^ 0, 
fi r admits exactly 2 W critical points (mod 2ir in each <j) n (x)) with energy bound (f21~f). 
The absolute minimizer of G(0, s) is obtained for 

!0, when sm ^P L ) > o, for every n = 1, . . . , N; 
7r, when sm( ^f^ < 0, for every n = 1, . . . , N, 
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By the Implicit Function Theorem argument in the previous paragraph, we conclude that 
for all r G (0, ri), Q r \s r is minimized by a unique s r = s(6i(r), . . . , Sj^(r)) with 



6 n (r) = 5* + 0{r), 



n — 1, 



AT, 



where 5* G {0,7r} is chosen as in (0). Finally, by Lemma |3.2| we conclude that for all 
r G (0,7*i) and sin(HpL) ^ the absolute minimizer of Q r in £ 9 is given by (f n ,(j) n ,A) = 
s r + rwi(s r ) + 0(r 2 ), with minimum energy 



inf SK = 2N P (l- r + 0(r 2 ). 



Hp 



When sin(HpL) = then Q r \s r i s degenerate at order r. Normally we should go to a 
higher order in the expansion to determine the stationary configurations at these exceptional 
values of H = rrnr/Lp, m = 1, 2, 3, . . ., but (as we will see later in the section) the existence of 
these degenerate points is explained by the exchange of stability of two families of solutions 
when vortices are nucleated into the sample from the boundary. We note that the treatment 
of the periodic problem in the subsequent paper [ [AlBeBr 00| | will require an expansion of 
the energy to order r 2 to resolve the degeneracy at any applied field H. 



4.2 Expanding the solutions to order r 

Because Q r is a smooth functional we may use the Implicit Function Theorem and the 



decomposition of Lemma [D] to obtain an expansion to arbitrary order in r of any critical 
point satisfying the energy estimate (Ell). Here we generate the expansion to order r, to get 
a better idea of the nature of the global minimizers. Take any such critical point of fl r \s r i 
with expansion as in (i) of Lemma pO, (f n , <f) n , A) = s + rwi(s) + 0(r 2 ), s = s(Si, . . . , Sn)- 



We deduce the equation satisfied by u>i(s) by implicit differentiation of the equation (|2 



-^P 1 - (\/n r (s + w(r,s))) 
dr 



P 1 f VT(s) + D 2 Q (s 
P 1 (VT(s) + D 2 Q (s) Ma)]). 



r=0 

dw 



(25) 

Since D 2 VLq{s) is an invertible map from W —>■ W this formula uniquely determines w\(s). 
Now we represent u>i(s) = (u n ,h v n,i, a x,u a z,i), m other words 



/„ = l + ru nil + 0(r 2 
A x = Hz + ra x i + 0(r 2 



A, 



a n + nHpx + v n> i + 0(r 2 ) 
ra z> i + 0(r 2 ), 
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and denote b(x, z) = curlai = d z a Xi i — d x a z> i. In terms of these coordinates ( pop takes the 
following simple form: 

^ ( cos(8 n + Hpx) + cos(5 n+ i + Hpx) — 2, 1 < n < N — 1, 



-— u n,i + 2m «,i = 2 i cos (^ + _ 1 

cos(5i + Hpx) — 1, 



n = N, 
n = 0, 



with <i(±L) = 0; 



1 
2 



' sin(5 n + i^pa;) — sin(5 n+ i + Hpx) — I n , 1 < n < N — 1, 
sin((5iv + -ffpx) — Jat, n = N, 

k — sin(5i + iJprr) — Jo, ft = 0, 



with boundary condition ^iL) — a X) i(±L, z n ) = and 



^ f /_z,[ s i n (<^i + -ffpx) — sin(5 n+ i + ifpx)] <ix, 1 < n < N — 1, 

n = 0; 



J_ L sin(5i + Hpx) dx 



and 



cti(x, z) 



dz 1 dx I 1 



At = b(x,z), ei 



where b(x, z) = b^ n \x) for z n -\ < z < z n , with 



dx ~ ~Y 



for n = 1, . . . , N. 



1 /" L 

sin(5 n + Hpx) / sin(5 n + Hpx) dx 

*2L J—L 



6 (n) (±L) = 0, 



(26) 



(27) 



Now assume that sin(HpL) ^ and consider the absolute minimizers, 5 n = 5* with 
5* G {0, 7r} chosen according to fl24]) . Except for an edge effect at the top and bottom layers 
(n = 0, iV) the gauge-invariant quantities are independent of n at order r. In particular, the 
magnetic field h, and Josephson current density are (at order r) completely independent 
of z, n: 



9 

K p . 



if j : = jcfnfn-i sin $ n , n _i = r— sin(5 + Hpx) + 0(r 2 ), 

h(x, z) = H + r— (cos(<5 + iJpi) - cos(<5 + Hpx)) + 0(r 2 ) 
2H 



(28) 
(29) 



where j c = is the critical Josephson current. The in-plane supercurrent vanishes at 
order r for all interior planes, 



J 



(n) 



0{r 2 



n = l,...,N-l, 



S. Alama, A.J. Berlinsky, & L. Bronsard 



20 



but the top and bottom of the sample carry current at order r: 

K 2 

jf) = + r— — (cos(5 + HpL) - cos(5 + Hpx)) + 0(r 2 ) (30) 
2Hp 

K 2 

jf) = + r— — (cos(5 + - cos(5 + HpL)) + 0(r 2 ) (31) 

Similarly, the order r expansion of modulus of the order parameter is n-independent on the 
interior planes, but is modified at the top and bottom at order r: 

f n (x) = 1 + r (~ + Acosh(^x) + B cos(6 + Hpx)\ + 0(r 2 ), n = 1, . . . , N- 1, (32) 



fo = In = 1 + ~ f ~ + A cosh(— x) + B cos(<5 + Hpx)) + 0(r 2 ) (33) 
2 \ 2 k / 



where 



n 3 Hp sin(5 + ZfpL) k 2 

^4 = — = 7— ~ i -D " 



V2(^p 2 + 2« a ) sinh(^) ' FV + 2k 2 ' 

Expansion to order r 2 in the solutions will show these quantities to be independent of z (or 
n) except for the top and bottom two planes and the top and bottom gaps. 
We note that the gauge-invariant phase difference, 

K 2 

$ n ,n-iO) = $ + Hpx + r— — (xHpcos(5 + HpL) - sin(5 + Hpx) + c n ) + 0(r 2 ), (34) 

2H Z 

contains constants c n which can only be determined by higher order expansion in r. In- 
deed, they represent the order r correction to the choice of 5 n when minimizing the finite 
dimensional problem r |5 r . 

We call this configuration vortex planes- see Figure 1. Unlike the Ginzburg-Landau 
case, the order parameter does not vanish at the "core" where the local field attains its 
maximal value and around which supercurrents circulate. The "vortices" are then the planes 
{x = const.} over which h(x) attains its relative maxima. These planes are nodes for the 
current, and flux per plane per cycle of the Josephson current is, to order r, the usual flux 
quantum {2n in our units.) 

Remark 4.2 (a) The distinction between the energies of the various lattice geometries at 
the lowest order term 0{r) in the expansion (|23|) is a surface term, in the sense that it scales 
like the length of the lateral edges of the sample 2Np, as opposed to the free energy itself 
which scales like the cross-sectional area of the bulk 2LNp. In other words, the vortex plane 
configuration is preferred for the effect it has on the surface currents on the left and right 
edges of the sample, with no regard to energy minimization in the interior. As we will see 
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in the periodic case in [[AlBeBr "OOfl the effect on the bulk will be observed in an order r 2 
term. For all other parameters fixed, eventually we can take r > small enough so that 
the 0(r) surface term dominates. But in a real superconductor the value of r is given, and 
hence when the width L is increased eventually the bulk 0(r 2 ) term will compete with the 
surface 0(r) term. This suggests that the interval of validity of the r expansion may not be 
uniform, but rather deteriorate with increasing sample width L. 

(b) Another indication that the radius of convergence of the expansion depends inversely on 
the width L is the presence of a linearly growing factor (secular term) in the order r term 



of the expansion fl34|) of $„, jn _i(x). 



4.3 Vortex plane nucleation 

The very precise desciption of minimizers for r ~ allows us to identify the transitions which 
the sample undergoes as new vortices are nucleated in an increasing applied field H. When 
< H < 7i /pL we observe the Meissner state for this problem: h attains its maximum at 
the boundary x = ±L, and decreases to a line of minima at x = 0. For our solutions the 
lower critical field Hex = Tr/pL: at this point the finite dimensional minimization problem 
G(0, s) — degenerates, and when H is increased slightly the minimizing configuration has 
phases d~ n = ti for all n = 1, . . . , N. Note that at this critical value of H a new node appears 
in the Josephson current at each endpoint. The switch in 5 n amounts to a flip in sign 
for the Josephson currents and for the variation of h from H, so the newly nucleated nodes 
correspond to minima of h, with the former minimum at x = becoming a local maximum. 
The same phenomenon will occur every time H crosses a value k = 2,3, . . . Each time, 
two new nodes for the Josephson current will be nucleated, and the change of 5 n from zero 
to 7i (or vice-versa) will exchange minima and maxima of h in the interior, resulting in the 
creation of exactly one new vortex plane. Since the minimum energy is given asymptotically 
by 

2 sm(HpL) 
Hp 

the magnetization M(H) := de/dH is discontinuous at each nucleation, indicating a first- 
order phase transition. 



e(H) := minQ r = Npr 2L 



+ 0(r 2 



5 The validity of the expansion 



In this section we prove Lemma |0| which justifies the finite dimensional reduction in a 
neighborhood of S. While parts (i)-(iii) follow easily from the Implicit Function Theorem, 
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the real interest is in the dependence of the interval of validity |r| < ro on the many param- 
eters in the model, especially the sample dimensions L and M = Np, k, and the applied 
field H. Here we give a lower bound for the "radius of convergence" in r as a function of 
these parameters. We show that the small r approximation is essentially independent of 
the number of planes N, but that its validity can be expected to decrease with increasing 
width L or increasing k. On the other hand, increasing the external field H enhances the 
approximations somewhat. This lower bound is consistent with experiments on the high-T c 
materials, where the vortex planes have not been observed, perhaps because the large values 
of k and the large size of typical samples (measured in terms of A a f>) reduce the radius of 
validity of the expansion below the appropriate value of r for such materials. For other types 
of layered superconductors with smaller values of k we are more likely to see configurations 
similar to the solutions produced in our r — > limit. 

Throughout this section, we concentrate on the parameters N, L, k, H, and assume < 
p < 1, k > 1, and L > 1. We also define a norm on the space W = [T5]- 1 : let w = 
(u n ,v n ,a) G W. Then we denote 

N rL 

IMP = \\[(Un,V n ,a)\\ 2 =pYj / (( U 'n) 2 + U 2 n + (v' n ) 2 + vty dx + 

n=0 J ~ L 

Our approach is to calculate a priori estimates on the solutions to the equation (|20|) to 
determine when that equation can degenerate. As a first step, we must recognize the system 
of differential equations satisfied by the solutions to (|20|) . For any U = (f n ,(f> n ,A) we write 
U = s + w(r, s), where s = s(Si, . . . , 5n) € S, and w(r, s) = (u n , v n , a) e W. We recall the 
effect that this decomposition has on some familiar quantities: 

/„ = 1 + u n , V n {x) = (j)' n - A x (x, z n ) = v' n - a x (x, z n ), 



\Va\ 2 + \a\ 



dx dz. 



a(x, z) = (d z £, -d x ti) , A£ = b(x,z), £\ 8 d = 0, h(x,z) — H + b(x,z), 



$„,n-l =4>n- 4>n-l ~ / A z (x, t) dt = HpX + S n + (fin^-!, 

[•Zn 

where ip n , n -i(x) = v n (x) - f„_i(x) - / a x (x,t) dt. 

J z„-i 

The equations for (u n ,v n ,a) will differ from the unconstrained Euler-Lagrange equations 
derived in section 2 because of Lagrange multipliers created by projection into the subspace 
W. 

To say w(r, s) solves P ± Z)f2 r (s + w(r, s)) = is equivalent to saying that Q r (s + w(r, s)) 
is stationary with respect to variations in W. Taking the variation with respect to just one 
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of the u n yields the exact same equations (0) as the unconstrained case, which we rewrite in 
terms of u„ : 



1 



K 



-u n + (u n + l)(u n + 2)u n + ^(v' n - a x (x, z n )) 



1 



u, 



K 



K 



>n - a x (x,z n )f 



(35) 



+ COS $„,„_! + (1 + u n +i) cos $„+!,„ - 2(1 + U n )) , n ^ Q, N; 



+ < 



((l + Ul)cOS$i j0 - (1 + Mo)) 



n = 0; 

71 = N, 



\ ((1 + Wjv-i) cos $jv,jv-1 - (1 + u N )) 

with boundary conditions u' n (±L) = 0. 

For v n , n = 0, . . . , N the integral constraint gives rise to the usual Lagrange multiplier 
without which the equations generally could not be integrated: 

{fnWn -a x (x,z n )) 



k 2 dx 



(36) 



|[/n/n-l SU1 $n,n-l ~~ fn+lfn Sm 3> n+ i )? J — C n 



f /at/at-i sin $jv,at-i _ Cat> 
-|/oAsin$i |Q + Co, 



n = 1, . 

n = N, 
n = 0, 



iV- 1; 



where 



Cm 



i /-£ | L/rJn-l Sin $ n>n _i - /„+!/„ Sm $ n+ i, n ] cftc, 71 = 1, 
f-L i/iv/jV-1 Sin $ N ,N-1 dx, 



AT-1; 



21/-L |/o/i sin $i )0 dx, 



1 r L r 
2L J-i 2- 

The equations for t> n also include a no-flux boundary condition, 



n = N, 
n = 0. 



/„K - «x(",2n)) 



a;=±L 



0, n = 0, . . . , N. 



Finally, we derive the Euler-Lagrange equations for a. As in our previous calculations it 
will be easier to deal with the associated magnetic field b(x, z) = curld(x, z). 

Lemma 5.1 If (f n , <f> n , A) = s + w(r, s) satisfies ([2d[), then b(x, z) = h(x, z) — H satisfies: 



b(x,z) = b (n \x), z n _ x <z<z n , n = l,...,N, 
6 W( X ) _ b {n+l \x) = p(v' n -a x (x,z n )), n = 1, . . . , N - 1, 
b {N \x) = p {(v' N - a x (x, z N ))fx) , b (1) (x) = -p [(v' - a x (x, z ))fo) 



db {n) 
dx 

where D„ 



rpn 
2 

1 rL 



/„/„_! sin $ n , n _ x (x) - D n , b^ n \±L) = 0, 



(37) 

(38) 
(39) 

(40) 



rpn 



2L J-l 2 



fnfn-l Sin $ n , n _i (x) dx, 
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Proof: First we derive the equation for a in the weak form. Any admissible variation for 
a in E can be represented in the form a = (d z rj, —d x rj) with rj G H 2 H Hq(D). (Indeed, rj 
solves Arj = curia in with T]\qd = 0: see Lemma [L|.) Taking the first variation of energy in 
the direction of a we obtain: 



— ( v 'n ~ a x(x, z n ))f%d z r}(x, z n ) dx 

n=u - L K 

+ £ / / \ijb(x,z)Ar] + ^-f n f n _ 1 sm^ njn _ 1 d x r]\ dzdx, (41) 

^J-LJzn-! Ik, 2 2 J 



for all rj G H 2 fl Hq(D). 

Next, we observe that there exist functions b n (x) G H l (—L,L) which solve the system 
(|38D-(|4"0D . Indeed, we define b n by integration, 

b n (x) = J l (^^-fn(x)f n ~i(x) sin ® n>n -i(x) - dx, n = l,...,N. 
Then b n (x) G H l ([—L, L]) and satisfies ([p. For n = 1, . . . , N - 1 we have: 

— (b n+1 (x)-b n (x)) = {fn+lfn Sin ® n +l,n- fnfn-1 Sm $ n ,„_i) + 

= -P^ (K - a ^(a;, *n))/*) , 

by (^). For the top and bottom gaps (n = 1, N) we obtain: 

d <tN/ \ rpn 2 Pf .- „ 
_6 (x) = —^-fNjN-ism9 N>N -i-D N 



N I ■ 



= V^j- ((v' N ~ a x (x,z N ))f 

^b°(x) = ^/o/isin$ li0 -Do 
ctx 2 

= -P^ (K - a x (x, z ))fo 

In each case we integrate (and use the boundary conditions b (-L) = = v' n (—L) 
a x (—L,z n )) to verify that b satisfies conditions (^) and 



Next we define b(x, z) = b (x) for z n -\ < z < z n , and show that b also solves (|4l|) . By 
an integration by parts in each summand, 

X/ / / —b(x,z)Ar]dzdx = 22 / / —b n (x)Ar]dzdx 

1 N ( [L ~ 
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/ — d x r](x, z) dz dx > 

-L J z n —i dx J 

— — - / 6 (x)c^77(x, ztv) — b(x)d z r](x, 0) c£r 
1 Ar ~ 1 ^ 

+ - E / f 6 "^) - d M%, Zn) dx 

~lL / ( ^/n/n-isin$„ iri _i(x) - £) n ) d x r]dzdx 

n= l J—L Jzn-l V 2 / 

= pH —{v' n -a x {x,z n ))fld z r]{x,z n )dx 

n=Q J-L K 

~Z2 / —fnfn-l^-^n,n-l{x)d x r]dzdx. 

for all T) G H 2 fl Hq(D). That is, 6 also solves the equation (^T|). Now consider the function 
<7(x, z) defined by g(x, z) = b(x, z) — b(x, z). Then g G L 2 (D) and solves 

^ g(x, z)Ar] dx dz = 

for all 77 G i? 2 fl Hq(D). Taking 77 to be the solution to the Dirichlet problem, Arj = g, 
v\dD = we arrive at the desired conclusion g = 0, and hence b(x,z) = b(x,z) and the 
lemma is established. 



As before, we observe that the gauge-invariant phase difference satisfies a useful equation 
(see (|10|)), which derives from Stokes' Theorem: 

-j-<Pn,n-l{x) = {v' n ~ d x (x, Z n )) - (v' n _ x - CL x (x, Z„_i)) + pb {n) (x) . (42) 

We now begin the proof of Lemma |3.1| . 
Step 1: Getting started. By the implicit function theorem, applied to the equation (|2H|), 
for every s G S, and for every fixed choice of parameters N, L, k, H, there exist constants 
p > 0, ri > such that (|2Up admits a unique solution w = w(r, s) with |r| < r x and \\w\\ < p. 
Since Q r (s + w(r, s)) is periodic in s we may compactify the problem by treating S as an 
iV-torus, and hence the constants r±,p may be chosen independently of s G S. Statements 
(i)-(iii) of Lemma |3.1| then follow from standard arguments involving the Implicit Function 
Theorem and the regularity of the functional Q r (see [|AmBa 98fl .) What remains to prove 
is the dependence of the value r x on the various parameters. 

We first observe that Lemma |3.1| (i) implies that the solutions w(r,s) satisfy a uniform 
(in s) estimate ||u>(r, s)\\ < Ci\r\, with constant G\ possibly depending on all parameters, 
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iV, L, k, H,p. From this preliminary argument we may, by reducing r% > if necessary, 
assume that 

\<M*)<\, n=0,...,N (43) 

for |r| < r\. 

Step 2: Uniform estimates in n. We now use the Euler-Lagrange equations fl3*5D-(£U) derived 
above to estimate the gauge-invariant quantities, in the supremum norm and individually in 
each plane or gap. 

First we make rough bounds, then iterate to improve the dependence on parameters. 
Since we can estimate the constants \C n \ < r and \D n \ < rn 2 p/2, we have 

||6 (rt) ||oc < crK 2 pL, (44) 

\\v' n - Oa;(-, 2V»)||cx3 < CTK 2 L. (45) 

(Here and in the following c will denote a pure constant, independent of all other parameters.) 

Now consider the equations for u n . Suppose first that max I ^«fc(x) = u n (xo) > 0. If 
n G {1, . . . , N — 1}, we observe that the extremal property of u n (xo) implies that 

u n -x(x ) cos (x ) COS 2u n {x ) < 0, 

while if n = or n = N the corresponding terms 

MAr_l(xo) COS$jv,AT-1 — Un(xq) < 0, Or Ul(xo) COS $1,0 — wo^o) < 0. 



Consequently, at the maximum point the right-hand side of (|3~5"D is non-positive. By the 
maximum principle we obtain u n (xo) < 0. (Note this argument holds in general, for any 
value of the parameters, for any solution of the Lawrence-Doniach system.) 

Next assume that min^ Uk(x) = u n (xo) < 0. By the same arguments as above at the 
point xq, the extremal u n satisfies an equation of the form 

1 1 



— zu'n(x ) + c(x )u n (x ) > — Av' n - a x (x , z n )) - 2r, 



with c(x) > |. In conclusion we may bound u n via 

2 

> Unix) > r sup(w^ - a x (x , z n )) 2 - Ar. (46) 

K z x 

Applying the simple estimate ( fI5D we obtain a preliminary bound on u n , 

> u n ix) > -or [tk 2 L 2 + l) . 

We use the calculus inequality 

2 5 
sup \f'ix)\<- sup + - sup \f"(x)\, for any < 5 < L, 

xe[-L,L] xe[-L,L] 1 xe[-L,L] 
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with 6 = 1, and obtain 

\\ u 'n\\oo < cr (rn 2 L 2 + 1 

To improve the above estimates we need to take into account the fact that the integrands 
and right-hand sides contain terms which oscillate rapidly when H is large. In particular, 
3?n,n-i — Hpx + (p ntn -i, and we expect that terms containing the sine or cosine of <& n ,n-i 
tend to zero in the weak sense as if — ► oo. Since these oscillatory terms are sources in 
the equations for the fields and currents, we should obtain stronger convergence to zero as 
H — ► oo, and hence sharper estimates for large H than we obtained in bluntly measuring 
these terms in the supremum norm. With this in mind, we require one further ingredient: a 
bound for tp' nn _ 1 using equation 



supK^a;)! <crn 2 L(l+p 2 ). (47) 

X 

We now integrate by parts in C n , D n , and in the integrals representing b^ n '(x) and 

bu P |6W(x)| < c^Cl + LlKIU + LII^IU (48) 

< c !^(i + rL V)(l + rL), 
Hp 



and similarly 



2 

TK 



sup\v' n (x) -a x (x,z n )\ < c-^(l + rLV)(l + rL). (49) 

Hp 



X 



Substituting into (HBI) we obtain a more refined estimate for u r 



2 



> u n {x) > -C u r [ 1 + rj^(l + rn 2 L 2 ) 2 (l + rL) 2 j , (50) 
with constant C u independent of N, k, L. 

Step 3: A lower bound on the linearization at r = 0. We define 

A = X(N,L,K,p) := ini{Q (V,V) : * = On, v n , a) G W, v n , a)\\ = 1}, 

with 

Qo(*,¥) = £> 2 O r=0 (/°,^,A°)[*,*] 

= PJ2 I (-^(/ i n) 2 + 2/i^ + ^(^-a :c (x,z„) 2 )) dx 

n=0 J ~ L ^ K K ' 

I rL i-Np 

+ ^7 / / (curl a) 2 dz dx. 

K J~L JO 
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Recall that, as an operator on the whole tangent space E, D 2 Qo(s) has an iV-dimensional 
kernel for any s G S. Here we consider its restriction to the orthogonal subpace W = TS^ = 
ker(D 2 r2 (s))- L . From the proof of Proposition |2J] (b) we may conclude that A > for any 
choice of N, L, K,p. Here we will obtain the more precise information on its dependence on 
these parameters. 

Note first that (trivially, since k > 1 by hypothesis,) 



N f L ( 1 A 1 N 



We require the following lemma to estimate the terms including the vector potential a 
appearing in the linearizations: 



Lemma 5.2 For any a G H 1 (D) with a ■ n — on dD, 



N 



srw ( Mia ^ (p + i)(JV + i) M l|2 

PL \\ a xV, z n)\\L 2 (D) — II a x II H^D) j 

n=0 iV 



AT 

E 



a z {z) dz^ 2 L 2( D } <p ll^a 2 z (x, z) dxdz. 



-\ J z n— 1 



(51) 



(52) 



Proof: Let x^(z) = (z — z n ~i)/p, n = 1, . . . , N. Then we apply the divergence theorem to 
the vector field (a x (x, z)xt,) 0) m the strip [— L, L] x (z n _i,2; n ) to obtain 



{a 2 (x, z n )) 2 dx 



< 



-L J z n -\ dz 

1 + - 

V 



(53) 



da T . 



dz 



dz dx 



for n = 
obtain 



' — L JZn-l 

N. Similarly, using Xn = ( z n+i — z)/p in the strip [— L, L] x (z n ,z n+ i) we 

2N 



[a x (x, z n )) 2 dx < 



1 

1 + - 
p 



-L 

for n = 0, ...,N- 1. 

Clearly there exists an index n such that 

2 



a 2 x + 



da x 



dz 



dz dx, 



(54) 



r Zn o 



a 2 + 



L J ZriQ — l 



da^ 



dz 



dzdx < — ,, 
~ N JJd 



a 2 + 



OOr, 



dz 



dz dx. 



For < n < uq — 1 we use estimate ([54]) for 1 1 ct-^ 1 1 § , and for no < n < N we use (|53|). In 
doing so we require the i7 1 -norm of a x in each interval exactly once, except for the strip 
.2 no _i < z < z no which appears twice. In this way we arrive at the desired estimate ([51]). 
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The estimate ([52]) is an elementary consequence of the Cauchy-Schwartz inequality 

We now estimate the second term of Q$: since 

2\v' n a x (x,z n )\ < (1 -e)(z/) 2 + ^—a 2 x (x, z n ), 

for any e > 0, we have 

(vL - a x (x,z n )) dx > e{u' n ) 2 dx - - a 2 x (x, z n ) dx. 

L J-L J-L 1 — e 

We now use ( |5T| ) and the elliptic estimate (|TT|)-(|I2"|) to estimate 

(curl a) 2 dxdz — p / al(x,z n )dx 

n =cr _L 1 ~~ e 

(if 4 L 2 N 2 p 2 Y 2 2e 

- ^2 l v 1 + ^4L2 + iw; (T^y 1 l|n|l//l 

>ari + ^r-^oai- 



2 V vr 2 J 1-e 



2^ ~ TTTTi N""* 1 ' 



uniformly in iV. Now choose e = |min{l, \ (l + -^L 2 ^ }• From the Poincare inequality 
(for H 1 functions on [-L, L] with vanishing mean) we may then obtain the lower bound 



j)V / [y' n — a x {x, z n )) 2 dx + // (curia) 2 dxdz 
n=0 J-L JJd 

N r L l f / 4 \~ 2 ) 

>epE j y n ) 2 dx + -mmy,(l + —L 2 ) j \\a\\ 2 m 

>e(l + ^L 2 ) V|:W^ + jniin|l,(l + ±L 2 ) '} \\a\\ 2 m 



-i min T' l 1+ ^ L ' 



We then obtain 

1 



■ AT 

IKIIhi + ll^ll^ 1 

. n=0 



A > — ^ min 1, ( 1 H — -L 2 ) \. (55) 



-3 

Note in particular that A is bounded away from zero uniformly in N, but the bound deteri- 



orates as either L — > oo or k — > oo. (See Remark 5^3 



Step 4: An upper bound on the linearization. We define Q r to be the quadratic form 
representing the second variation of energy around the solution s + w(r, s) G <S r , with respect 
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to variations $ = (fi n , u n , a) in the subspace W: 

tt r (s + w(r,s) +tV) 



2 dt 2 



t=o 



P E / + -K - + (3/„ 2 - l)/i 2 

4 1 

H — ?( v 'n ~ ax)fnHnW n - a x (x, z n )) + —= f 2 (u' n - a x (x, z n )) 2 \ dx 



+ o P E / 1 (fn^n-l + fn-lfin) sin $n,n-lK ~ f n -l ~ / «2 ^] 
fnfn— 

1 COS $ n , n -lK - ~ / 

J Zn—l 

+/4 + ~ 2 Vn^n-l COS $ n , n _i } rfx 

+ — // (curia) 2 dxdz. 
k JJb 

As long as the linearization Q r remains non-singular (as an operator on W), we may 
invoke the Implicit Function Theorem to ensure that the solutions to equation (EDI) in the 



form s + w(r, s) determined by Lemma 3J. remain valid. Hence, the critical value r$ which 
determines the radius of convergence of the expansion is bounded below by the smallest 
value of r for which Q r admits zero as an eigenvalue. If Q r degenerates at r*, then there 
exists a test vector \]> e W, ||\&|| = 1, such that Q rt (ty, ^f) = 0. Using the lower bound on 
Qo, we have: 

A < Q (*, *) = Q r , (vl/, *) + (Q - Q r J (% *) = (Q - Q r J (tf , (56) 

An estimate of the difference between the linearization at r — r* and at r = in terms of 
r* (and the norm of \1/) will then yield a lower bound on the critical value tq > r*. This 
estimate follows from the sup-norm bounds on each u n , (v' n — a x (-, z n )) obtained in ( EUJ ) and 

(H): 

N r L r 1 

IQ r .(*,*)-Q (*,*)l < -K-^)V 2 + 3|/ 2 -i|/i 2 

n=0 -/-£ U 

4 / / 



K 2 



AT ,L 



+ — pJ2 / i(l)"n-l| + - / " 2 ^z| 

+ - 2|AW*n-i| + K - y n-i - / « 2 dz} 2 \ dx 



< cr* [1 + K(l + nK 2 K)} || ^ 11 2 



S. Alama, A.J. Berlinsky, & L. Bronsard 



31 



where 

K:= _(l + r *LV)(l+r*L). 

Inserting the above estimate and (|5"5j ) into ( pUD produces a lower bound on the first point of 
degeneracy r*, 

r* [1 + #(1 + r*K 2 K)] > A > ^ min jl, (l + ^L 2 ) " j . (57) 

We observe that the left-hand term is monotone increasing in r*, and therefore implies a 
lower bound of the form 

ro> > R(L, k, H), 

where the constant R can be chosen uniformly in N, but decreases with increasing L and k 
and increases with increasing H. 

Step 5: Conclusion. The above estimates have all been based on the (unsatisfactory) initial 
hypothesis \r\ < ri(N,L,K,H) chosen small enough to ensure (f43[). However, we observe 
that the only role of this hypothesis is to obtain fl43|), and the estimates obtained thereafter 
are valid whenever (^) holds. Since ||w n ||oo is a continuous function of r, the estimates 
obtained for u n ,v n ,a and r* above persist as we increase r, either until (pE3| ) is first violated 
or until we reach r = r*. From fl50|), if condition (^3|) is violated then r satisfies 

C u r[l + rK 2 K 2 } > 1 

a condition on r which also defines a lower bound which is independent of N, decreasing 
with increasing L, k, and increasing as H increases (just as for R(L, k, H) above.) In either 
case, the perturbation argument remains valid for an interval r e [0, tq) with ro bounded 
below by a constant which is independent of N, with the dependences on L, k, H claimed in 
the statement of (iv). This concludes the proof of Theorem |3.1| . 
<> 



Remark 5.3 We note that the lower bound on r* obtained above improves with smaller 
L, k and larger H. Of course this is only a one-sided bound, and we cannot be sure that 
the solutions obtained by the degenerate perturbation method really do cease to exist if L 
or k are too large. Indeed, solutions with a similar form may persist beyond the range of 



validity of Lemma [3.1|, and the estimates which gave the lower bound (|57|) are certainly not 
sharp. However, we can show that the smallest eigenvalue A (AT, L, k) of the linearization at 
r = does tend to zero as either L, k —>■ oo. From the Implicit Function Theorem we may 
then infer that w(r, s) grows rapidly near r = 0, which in turn suggests that the interval of 
validity is diminished as L, k increase. 
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Lemma 5.4 For any fixed N, 

as L, k — > oo. 

Proof: As before, set M = 7Vp. Let £ G H 2 \~\Hl(D) be the solution to the Dirichlet problem 

A£ = l, e|aD = 0, 
and define a = (d z £, —d x £). Note that £ is even in x, and 

— (I (cm\ a) 2 dxdz 

K J J D 



2LM 



We next choose v n with 



d f L 



J v n {x) dx = 0. 



dx 

Since ^(i, z) is even in x we have 

rx 

v n {x) = / a x (x', z n ) dx', 
Jo 

and z/„ is odd. Applying Stokes' Theorem to the rectangle [—L,L] x [z n -i,z n ], 

rz n rx ^ 

2px — I / curl a dxdz 

JZn-l J-X 

/x rz„ 
(a x (x, z n ) - a x (x, z n _i)) dx + (a z (-x, z) - a z (x, z)) dz 
'X J Z„-l 

= 2 {v n {x) - v n -i(x) - J a z (x, z) dzj . 

In particular, 

N N 

2^2\Wn\\l > Yl [\\ U n\\l + \Wn-l\\l 

n=0 n=l 

9 n=l J ~ L V^-l J 

> — -^—L 3 —p II [a z (x, z)} 2 dz dx. 
9 JJd 

Taking ji n = 0, \I/ = (fj, n , v n , a), we have 

2 \ \ ~> II ||2 i || -»||2 \ 3 
E>P2^ IWnh + IM| 2 > —^~ L ' 
n=0 9 
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Figure 1: Vortex planes, for a sample with a finite number of superconducting planes. The dotted 

horizontal lines represent the superconducting planes, and horizontal arrows indicate the in-plane 

in) 

currents jx , which vanish to order r in the interior of the sample. The vertical arrows indicate 

in) 

the Josephson currents j\ in each gap. Except for an edge effect near the top and bottom of 
the sample, the magnetic Geld h(x, z) and Josephson currents j z are z-independent. The vortices 
correspond to local maxima of h, and lie along planes x =constant, indicated here by the dark 
bands. 



assuming (as usual) that < p < 1. We may then conclude that: 

Qo(*,¥) 



X(N,L,k) < 



1*111 



< ^ Id (curl oT) 2 dxdz 



9 

- 2k 2 P 2 L 2 ' 
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